1. Introduction
Global navigational satellite systems (GNSS) include Global Positioning System (GPS) (United States), GLONASS (Russia), Galileo (Europe) and COMPASS (China) (systems in use or in development). A GNSS typically uses a plurality of satellites orbiting the earth. The plurality of satellites forms a constellation of satellites. A GNSS receiver detects a code modulated on an electromagnetic signal broadcasted by a satellite. The code is also called a ranging code. Code detection includes comparing the bit sequence modulated on the broadcasted signal with a receiver-side version of the code to be detected. Based on the detection of the time of arrival of the code for each of a series of the satellites, the GNSS receiver estimates its position. Positioning includes geolocation, i.e. the positioning on the surface of the Earth.
An overview of GPS, GLONASS and Galileo is provided in sections 2.1.1, 2.1.2 and 2.1.3 of Sandra Verhagen, The GNSS integer ambiguities: estimation and validation, Delft University of Technology, 2004, ISBN 90-804147-4-3 (herein referred as “[1]”) (which was also published in Publications on Geodesy 58, Delft, 2005, ISBN-13: 978 90 6132 290 0, ISBN-10: 90 6132 290).
Positioning using GNSS signal codes provides a limited accuracy, notably due to the distortion the code is subject to upon transmission through the atmosphere. For instance, the GPS includes the transmission of a coarse/acquisition (C/A) code at 1575.45 MHz, the so-called L1 frequency. This code is freely available to the public, in comparison to the Precise (P) code, which is reserved for military applications. The accuracy of code-based positioning using the GPS C/A code is approximately 15 meters, when taking into account both the electronic uncertainty associated with the detection of the C/A code (electronic detection of the time of arrival of the pseudorandom code) and other errors including those caused by ionospheric and tropospheric effects, ephemeris errors, satellite clock errors and multipath propagation.
An alternative to positioning based on the detection of a code is positioning based on carrier phase measurements. In this alternative approach, the carrier phase of the GNSS signal transmitted from the satellite is detected, not the code modulated on the signal transmitted from the satellite.
The approach based on carrier phase measurements has the potential to provide much greater position precision, i.e. up to centimeter-level or even millimeter-level precision, compared to the code-based approach. The reason may be intuitively understood as follows. The code, such as the GPS C/A code on the L1 band, is much longer than one cycle of the carrier on which the code is modulated. The position resolution may therefore be viewed as greater for carrier phase detection than for code detection.
However, in the process of estimating the position based on carrier phase measurements, the carrier phases are ambiguous by an unknown number of cycles (this is for instance explained in [1], section 1.1, second paragraph). The phase of a received signal can be determined, but the cycle cannot be directly determined in an unambiguous manner. This is the so-called “integer ambiguity problem”, “integer ambiguity resolution problem” or “phase ambiguity resolution problem”.
GNSS observation equations for code observations and for phase observations are for instance provided in [1], respectively sections 2.2.1 and 2.2.2. An introduction to the GNSS integer resolution problem is provided in [1], section 3. The idea of using carrier phase data for GNSS positioning was however already introduced in 1984 in Remondi, Using the Global Positioning System (GPS) Phase Observable for Relative Geodesy: Modeling Processing and Results, Center for Space Research, The University of Texas at Austin, May, 1984 (herein referred as “[2]”).
The basic principles of the GNSS integer resolution problem will be now explained with reference to FIGS. 1 to 4. Further explanations are then provided, with mathematical support and explanations of the further factors generally involved in implementing an integer resolution system for precise position estimation.
2. Basic Principles of the GNSS Integer Resolution Problem
FIG. 1 schematically illustrates a GNSS with only two satellites and one receiver. A mobile receiver is also called a rover, a stationary receiver may be called a base station (useful for differential processing, e.g., DGPS) or a reference station (useful for network processing). While an actual GNSS involves more than two satellites, in order to illustrate the basic principles of GNSS carrier-phase measurements, only two satellites are represented on FIG. 1. Each satellite broadcasts a signal (illustrated by the arrows originating from the satellites and going in the direction of the receiver). A portion of the carrier of each of the signals is schematically represented in the vicinity of the receiver. The portion of the carrier is represented in FIG. 1 in a non-modulated form for the sake of clarity. In reality, a code is modulated, e.g. bi-phase shift key (BPSK)-modulated, on the caner.
At one point in time, the receiver can measure the phase of the carrier of the received signal. The receiver can also track the carrier phase, and lock onto it to track the additional cycles of the carrier due to the changing distance between the receiver and the satellite. The aspect of tracking the phase in time is however disregarded for now, for the sake of clarity, but will be explained later.
While the receiver can measure the phase of the carrier of the received signal, the number of cycles between the satellite and the receiver is unknown. The distance corresponding for instance to one cycle of a GPS L1 frequency carrier (1575.45 MHz) is about 19 centimeters (taking into account the propagation speed of the radio signal, i.e. the speed of light). The distance of 19 centimeters corresponds to the wavelength of the carrier. In other words, the carrier phase can be measured starting at one point in time, but the integer number of cycles from the satellite to the receiver, e.g. at the starting point in time, is unknown. Also, the difference between the number of cycles from the first satellite to the receiver and the number of cycles from the second satellite to the receiver is unknown.
The amplitude peaks of some successive cycles of a GNSS signal carrier in FIG. 1 are illustrated by the indications “−1”, “0”, “+1”, “+2”. This intuitively illustrates the integer ambiguity as to the number of cycles. The peak referred to by the indication “0” corresponds to a particular number of cycles from the receiver to the satellite. The peak referred to by the indication “+1” corresponds to one more cycle, the peak referred to by the indication “−1” corresponds to one less cycle, and so on. The estimation of the position of the receiver depends on the knowledge of the correct number of cycles from the satellite to the receiver.
A given number of cycles from the first satellite to the receiver (e.g. corresponding to the indication “+1” on the left-hand side of FIG. 1) and a given number of cycles from the second satellite to the receiver (e.g. corresponding to the indication “−1” on the right-hand side of FIG. 1) leads to one estimated position. In other words, it can be seen that the integer couple (+1, −1) leads to one estimated position in the two-dimensional illustration. Another combination of number of cycles to each of the satellites leads to another estimated position. For instance, the combination (0, +1) leads to another estimated position in the two-dimensional illustration. This will be better understood with reference to FIGS. 2a to 4.
FIG. 2a discloses a two-dimensional uncertainty region wherein the receiver is known to be, or at least is known to be very likely to be. The knowledge that the receiver is within this uncertainty circle may for instance be based on a coarse code-based position estimation. While an actual GNSS implementation involves more than a two-dimensional uncertainty region, the basic principles of carrier-phase measurements may be understood by referring to a two-dimensional space of unknowns only. It is known or assumed that the position of the receiver is within the uncertainty region limited by the represented circle, but the position within the circle is unknown yet.
FIG. 2b discloses the two-dimensional uncertainty region of FIG. 2a with, in addition, a dashed line with an arrow. The dashed line represents the direction towards one GNSS satellite, called Satellite 1. Satellite 1 broadcasts a signal containing a carrier. The parallel lines represented in the uncertainty region are the lines wherein the receiver may be located, when assuming that only the carrier of Satellite 1 is used for position estimation. Each line corresponds to a particular wave front of the signal broadcasted by Satellite 1. Adjacent wave fronts are separated by one carrier wavelength. Depending on the correct number of cycles to Satellite 1, the receiver may be viewed as being located on one of these lines. Due to the integer ambiguity, the correct position line is unknown.
FIG. 2c represents the same elements as FIG. 2b, i.e. the uncertainty region and the carrier wave fronts of the signal broadcasted by Satellite 1. In addition, FIG. 2c includes a second dashed line with an arrow schematically indicating the direction of a second GNSS satellite, i.e. Satellite 2. The lines which are perpendicular to the direction of Satellite 2 represent the carrier wave fronts of the signal broadcasted by Satellite 2. Depending on the correct number of cycles to Satellite 2, the receiver may be viewed as being located on one of these lines.
The carrier of the signal broadcasted by Satellite 2 may be used for position estimation in addition to the position information derived from the carrier associated with Satellite 1. The uncertainty can therefore be reduced. Instead of being able to only assume that the position of the receiver is somewhere within the circle (of FIG. 2a), or to only assume that the position of the receiver is on one of the straight lines represented on FIG. 2b, it can now be assumed that the position of the receiver is on one of the lines corresponding to the wave fronts of Satellite 1 and simultaneously on one of the parallel lines corresponding to Satellite 2. In other words, by using carrier phase measurements from the two satellites, it can be assumed that the position of the receiver is or is very likely to be on the intersection between a line associated with a wave front of Satellite 1 and a line associated with a wave front of Satellite 2.
It can be seen that there are 28 line intersections within the exemplary uncertainty region of FIG. 2c. These intersections correspond to 28 possible positions for the receiver. The position estimation problem within the uncertainty region is therefore reduced to the problem of finding on which one of these intersections the receiver is located.
More than two satellites are available to estimate the receiver position. The carrier of a signal from a further satellite can therefore be used in an attempt to resolve the integer ambiguity.
FIG. 2d represents the same elements as FIG. 2c, i.e. the uncertainty region and the carrier wave fronts of the signals broadcasted by Satellite 1 and by Satellite 2. In addition, FIG. 2d includes a third dashed line with an arrow indicating the direction of a third GNSS satellite, i.e. Satellite 3. Accordingly, a third set of carrier wave fronts, associated with Satellite 3, are represented.
The third set of wave fronts added on FIG. 2d and superposed on the pattern with the exemplary 28 intersections of FIG. 2c provides additional information to assist in resolving the integer ambiguity. How well a wave front associated with Satellite 3 fits with an intersection between the wave fronts associated with Satellite 1 and Satellite 2 provides an indication as to the probability that one particular intersection is the correct position. For instance, the triplet (0, 0, 0) can be intuitively regarded as providing a highly likely combination of wave fronts. The receiver is likely to be located at the intersection corresponding to the triplets (0, 0, 0) of wave fronts. This is however not the only possible intersection, also called a node.
In order to assign to each intersection of wave fronts associated with Satellite 1 and Satellite 2 (as illustrated on FIG. 2c) a probability of being the correct position, each intersection may be considered one by one. FIG. 2e illustrates a particular intersection, or node, being considered, or searched. The particular intersection, or search node location (using the label of FIG. 2e), is (+2, 0), i.e. wave front “+2” associated with Satellite 1 with wave front “0” associated with Satellite 2. How far the search node is from the closest wave front associated with Satellite 3 may intuitively be viewed as providing an indication of the probability that the search node corresponds to the correct position.
The wave front “+2” associated with Satellite 3 is the closest wave front from the search node (+2, 0). The wave front “+2” associated with Satellite 3 is however relatively far from the search node, i.e. the superposition match is not good. It results that the triplet (+2, 0, +2) is relatively unlikely to correspond to the correct position.
Each combination of a wave front associated with Satellite 1 and a wave front associated with Satellite 2 may be considered and a probability may be assigned to each one of these integer combinations, also called “ambiguities” or “integer ambiguities”. This is shown in FIG. 3. The vertical arrow which is associated with each intersection, or search node, gives an indication of the probability that the search node corresponds to the correct position. Note that the probabilities shown on FIG. 3 are exemplary and do not necessarily correspond to the situation illustrated in FIGS. 2d and 2e. 
FIG. 3 shows that the couple (0, 0) has a predominant probability within the uncertainty region. It is highly likely that this couple of integers corresponds to the correct integer solution from which the receiver position can be most precisely derived.
FIG. 4 shows another exemplary probability map wherein the probability associated with each intersection of wave fronts is indicated by a vertical arrow. However, in FIG. 4, two search nodes, corresponding to the couples (0, 0) and (+1, −1), have an almost equal probability. In this situation, it is highly likely that one of these two integer couples corresponds to the correct integer solution from which the receiver position can be most precisely derived.
The GNSS integer resolution problem is therefore the problem consisting in estimating with the highest possible confidence the correct set of integer values for the ambiguities, so as to provide high precision GNSS positioning. It may include scanning through possible integer combinations or search nodes (as illustrated in FIG. 2e) using statistical measures to assess which combination is the correct one.
3. Mathematical Formulation and Further Considerations
3.1 Observation Equations
The linear GNSS observation equations are provided byy=Aa+Bb+e  (1)wherein                y is the GPS observation vector of order m (the vector of observables),        a and b are the unknown parameter vectors of dimension n and p respectively,        A and B represent the design matrices derived from the linearized observation equations of the GNSS model used, and        e is the noise vector (the residuals).        
The entries of the parameter vector a are unknown integer carrier phase ambiguities, which are expressed in units of cycles, i.e. aεZn. The remaining parameters are the so-called baseline parameters, i.e. bεRpincluding for instance atmospheric delays (see [1], sections 2.2 and 3.1). When the number of observations increases, the noise vector becomes to behave in accordance with a normal distribution.
The observation vector y may comprise both phase and code observations, on as many frequencies as available, and accumulated over all observation epochs. The problem described in above section “[2. Basic principles of the GNSS integer resolution problem]” is therefore generalized with many observations, based on one or more frequencies from a plurality of satellites for determining the position of a receiver.
The observation equations may take into account many observation data and many types of observation data. For instance, differential data obtained from a differential GPS (DGPS) system. DGPS uses one or more reference stations whose positions are precisely known. This enables to compute the effect of the ionosphere, satellite clock drifts, satellite ephemeris error at one moment in time. DGPS techniques were originally developed to compensate for the intentional clock error added when the so-called “Selective Availability (SA)” feature of the GPS was turned on. The SA feature has been permanently turned off in 2000, but DGPS techniques are still useful to compensate for effects which cause a reasonably constant delay during a certain period of time. The information can therefore be provided to the receiver, which can take it into account for improving the position estimation precision. In other words, the DGPS reference stations transmit differential corrections to be used by the GNSS receiver to improve the position estimation. These corrections may be integrated as observations into the observation equations without affecting the general principles behind the GNSS integer ambiguity resolution problem.
The system of equations (1) and the problem associated with its resolution are particular for the following reason. The equations involve real number unknowns (baseline parameters b) and integer number unknowns (integer ambiguities a). Techniques are not well-established to handle the resolution of a system of equations where some of the unknowns are known to have integer values.
In the above observation system (1), the number of unknowns is not equal to the number of observations. There are more observations than the number of unknowns, so that the system is an overdetermined system of equations. Statistical properties can be associated to each possible solution of the overdetermined system. The resolution of the system consists in finding the most likely solution to the system, including values for the integer number unknowns and real number unknowns.
One specific characteristic of carrier phase measurement based GNSS models is that, if an integer solution can be identified with a high confidence (i.e. with a probability of correct solution close to one), this provides a very precise solution (up to centimeter-level precision). However, integer ambiguity resolution carries the risk that the carrier phase ambiguities on one or more satellites and frequency bands are incorrectly fixed. Because GNSS GPS signals have a carrier of about 20 centimeters, if a carrier phase ambiguity is set to a wrong integer, this may result in position errors of decimeters or more.
3.2 Generalization in Time, and Float Solution
In above section “[2. Basic principles of the GNSS integer resolution problem]”, only one snapshot in time was considered. However, the precision may be improved over time, by using a series of successive observations.
The phase of the carrier broadcasted from one satellite differs from one observation at a first point in time to another observation at a second point in time. However the carrier phases may be tracked, so that the carrier phase ambiguities themselves do not change. In other words, the receiver can lock onto the phase of one carrier.
A usual technique to solve the system of equations (1) is to first treat the unknowns a as real numbers, i.e. floating-point numbers (even though the unknowns of the vector a are known to be integer numbers). This approach is advantageous because well-known techniques such as Kalman filtering and least squares resolution may be used to derive a float solution, i.e. a set of real numbers for the unknown carrier phase ambiguities. This approach is advantageous even though it does not take into account the fact that the unknown carrier phase ambiguities, i.e. the entries of the parameter vector a in the equation system (1) above, are integer values.
The approach is also advantageous because it converges to the correct integer solution, provided that sufficient time is allowed. By adding snapshots (sets of observations at respective epochs), and by keeping observations associated with each snapshot to improve the solution (e.g. least-squares solution) a converging solution emerges.
However, a practical disadvantage of the float solution is that, while it converges, it takes time to do so. Typically, float solutions generally takes tens of minutes to converge sufficiently for centimeter-level work. While the Kalman filter may be used to obtain a float solution, and the float solution converges to a good estimate of the carrier phase ambiguities, it takes a long time to do so. There is a need to speed up this process.
The fact that the unknowns a are integer numbers may be used as constraints, for deriving combinations of integer values from the float solution, even before the float solution has converged sufficiently to the correct integer solution. This will now be explained.
3.3 Fixing the Float Solution to Integer Values: Fixed Solution
Known approaches exist to fix the ambiguities to integer values. For instance, the float solution is projected or mapped to an integer solution. See sections 2.1 and 2.2 of Teunissen, P. J. G. (2003), GNSS Best Integer Equivariant estimation, presented at IUGG2003, session G04, Sapporo, Japan (herein referred as “[3]”).
This enables to obtain a fixed solution quicker. The reason behind fixing the float solution to an integer solution is to use the knowledge that the ambiguities in fact must be integer so that, if the correct integers are selected based on the float values of the ambiguities, the estimation of the remaining unknowns (i.e. the unknowns of vector b in equation system (1)) is improved and accelerated. The ambiguities are therefore fixed to the most probable integer solution, to reduce the number of unknowns and to thereby increase the over-determination of the system of equations.
Fixing includes setting each ambiguity to an integer value and using the combination of integer values to reduce the number of unknowns into the system of equations. However, while fixing the ambiguities speeds up the convergence process, this comes at the risk of incorrectly fixing the ambiguities, which can lead to the convergence towards a wrong solution.
It is noted that the ambiguities indicate the unknown number of cycles of the carrier between the receiver and the satellite at a particular instance in time, e.g. when initializing the system or at any other starting point in time, and thus are fixed values. The filter thus is used to estimate these fixed values based on the observations made. More precisely, by collecting more and more observations over time, the state vector of the filter, comprising inter alia the ambiguities, gradually converges to stable integer values representing the ambiguities.
To account for the changing distance from the satellite to the receiver over time usually a phase locked loop is employed in the receiver to track the carrier signal to determine the additional number of cycles to be added or deducted from the initial value that is to be estimated by the filter for the ambiguities at the starting point in time.
The Best Integer Equivariant (BIE) estimator (see section 4 “Best Integer Equivariant Estimation” of [3], or chapter 4 of [1]) is an example of estimator in which the integer-nature of the carrier phase ambiguities is used without explicitly enforcing a single correct integer ambiguity combination. The BIE approach uses a weighted average of integer ambiguity combinations to produce a solution which has a precision that is always better than or as good as the precision of its float and fixed counterparts (see [1], section 4.1, page 69).
An important computational aspect of the BIE approach is that the float ambiguity solution should be transformed into a more orthogonal space to accelerate the process used to generate integer ambiguity combinations used in the weighted summation (see [1], section 4.2.2, page 71). This transformation into a more orthogonal search space is referred to as “Z-transformation” and is explained in more details for instance in [1], pages 33-36; section 3.1.4 or in Teunissen, The least-squares ambiguity decorrelation adjustment: a method for fast GPS integer ambiguity estimation, Journal of Geodesy, 70: 65-82. A two-dimensional example of Z-transformation of an ambiguity search space is shown in [1], page 34, FIG. 3.6, in page 34. As mentioned in [1], page 34, lines 3-6: “Due to the high correlation between the individual ambiguities, the search space in the case of GNSS is extremely elongated, so that the search for the integer solution may take very long. Therefore, the search space is first transformed to a more spherical shape by means of a decorrelation of the original float ambiguities.”
The Z-transformation, also called herein simply Z-transform, is different from (and should not be confused with) a Fourier transform in the frequency domain. The Z-transform is also different from (and should not be confused with) the mapping of the float values from a Rn real number space to a Zn integer number space.
4. Weighting Over the Ambiguities
Methods have also been disclosed to estimate the receiver position based on all possible integer combinations.
Betti B., Crespi M., Sanso F., A geometric illustration of ambiguity resolution in GPS theory and a Bayesian approach, Manuscripta Geodaetica (1993) 18:317-330 (herein referred as “[4]”) discloses a method wherein the ambiguities do not need to be resolved, but rather the method involves “sum[ming] over all possible ambiguities with proper weights directly derived from the likelihood function” (page 326, left-hand column, Remark 4.2).
In [4], it is also suggested to restrict the averaging to some of the combinations of integer values: “It has to be underlined that on the practical implementation of . . . in reality we have extended the summation not over the whole grid of ambiguities but just to the closer knots as the function . . . drops very quickly to zero when β attains large values” (page 327, right-hand column, lines 34-39).
5. Problem to be Solved
There is a need for improving the implementation of positioning systems based on GNSS carrier phase measurements, to obtain a precise estimation of the receiver position in quick, stable and user-friendly manner.